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ABSTRACT 

We describe a three-dimensional simulation of a \Mq solar-type star approaching a lO^M© black hole on 
a parabolic orbit with a pericenter distance well within the tidal radius. While falling towards the black hole, 
the star is not only stretched along the orbital direction but even more severely compressed at right angles to 
the orbit. The overbearing degree of compression achieved shortly after pericenter leads to the production of 
strong shocks which largely homogenize the temperature profile of the star, resulting in surface temperatures 
comparable to the initial temperature of the star's core. This phenomenon, which precedes the fallback ac- 
cretion phase, gives rise to a unique double-peaked X-ray signature that, if detected, may be one of the few 
observable diagnostics of how stars behave under the influence of strong gravitational fields. If ~ lO^M© black 
holes were prevalent in small or even dwarf galaxies, the nearest of such flares may be detectable by EXIST 
from no further away than the Virgo Cluster. 

Subject headings: X-rays; bursts — shock waves — black hole physics — radiation mechanisms: thermal 
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1. INTRODUCTION 

Each star in the vicinity of a supermassive black hole 
(SMBH) traces out a complicated orbit under the combined 
influence of all the other stars and the black hole itself. There 
is a chance that encounters with other stars shift a star onto a 
nearly radial loss cone orbit which brings it very close to the 
black hole. Exactly how close a star can approach a black hole 
without suffering distortions is defined by the tidal radius 

rt ~ 7 X IQ^^m]!^ {R*/Rq) (M^/Mq)-'/^ cm, (1) 

where M(, denotes the mass of the black hole in units of 
IO^Mq. While the tidal radius is defined by average stellar 
properties, a solar-type star passing within this distance would 
likely be disrupted in a single flyby. 

It is a complicated (although tractable) problem of stellar 
dynamics to calculate how frequently a star enters this zone 
of vulnerability ( Frank & Rees 1976;). For galaxies with steep 
density cusps, such disru ptions would take place about once 
every few thousand years (iMagorrian & Tremaind 19991) . The 
exact rate depends on the statistics of the stellar orbits and 
particularly on how quickly the near-radial loss cone orbits 
are replenished. When a star is disrupted, the sudden release 
of gravitational binding energy is bound to produce a burst of 
radiation. The flares resulting from a disrupted star could be 
the clearest diagnostic of a black hole's presence. 

When a rapidly changing tidal force starts to compete 
with the star's self-gra vity, the ma t erial of the star responds 
in a co mplicated way ( Reeslll988l: ICarter & Lumineti [1981: 
iBicknel l & Gingold 1983). During a close passage, the star is 
stretched along the orbital direction, squeezed at a right angle 
to the orbit, and strongly compressed in the direction perpen- 
dicular to the orbital plane. This phenomenon poses a difficult 
challenge to computer simulations — three-dimensional gas- 
dynamical calculations have so far addressed the fate of the 
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bulk of the matter, but key questions relating to the details of 
extreme tidal compression have yet to be answered. In par- 
ticular, very high spatial resolution is needed to model stars 
passing well within the tidal radius, for which extreme com- 
pression is halted by a shock which rebounds and eventually 
breaks out of the stellar surface. In this paper, we present 
the most highly resolved three-dimensional simulation of the 
tidal disruption of a solar-type star to dat e, with voxels ^10^ 
times smaller than previous calculations (iLodato et al.ll2009t 
iKob avashi et al. 2004). We use the results of our simulation 
to construct a model for the breakout of Shockwaves from a 
tidally disrupted star, and we apply this model to determine 
the detectability of these tidal shock breakouts (TSBs) in the 
local Universe. Our results provide new insights into what 
happens when stars are strongly shocked as a result of ex- 
treme compression and the characteristic properties of such 
events. 

The structure of the present article is as follows. A descrip- 
tion of the numerical methods and the initial models are sum- 
marized in ^ In ^ we present a chronological overview of 
the events in the life of a tidally disrupted star An analyti- 
cal treatment of shock breakout in tidally disrupted stars fol- 
lows in ^ while estimates for the luminosity and detectabil- 
ity of various TSB events are subsequently presented in ^ 
We identify limitations and summarize our results in ^ For 
convenience, a glossary of symbols (Table [T]) is included as 
an appendix. 

2. NUMERICAL METHOD AND INITIAL MODEL 



Our simulation is carried out using FLASH jFrvxeU et al.l 
|2003), an adaptive mesh code that has been used to treat a 
wide variety of gas dynamics problems. The black hole is ini- 
tialized as a point mass with Mh = lO^M©, while the progeni- 
tor star is initialized as a 1 Mq, F = 3 polytrope (P cx p'"'"''''") 
with /?*=7xl0"^ cm and central density of 76 g cm""*. Dur- 
ing the simulation, the gas obeys a 7ad = 5/3 equation of state. 
To ensure that the hydrostatic equilibrium of this configura- 
tion is maintained, we ran a control simulation without the 
presence of the black hole for lO'* s, roughly the duration of 
the encounter We found good stability — density does not 
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(c) r = s 



Fig. 1. — Events in the life of a tidally-disrupted solar-type star. The star travels from left to right in the black hole's frame, with f = s when r = ip. the 
pericenter passage distance. The event horizon, which lies a distance from the black hole, is shown as a thick dashed line. The orange coloring represents 
logjQ p, and the light blue line shows the path of the star's center of mass. The encounter is not drawn to scale. 



change by more than a part in 100 for the inner 99.5% of the 
star's mass. 

The center of the computational domain is fixed to the star's 
center of mass and is 8 x 10" cm on a side. We subtract out 
the gravitational force the black hole exerts on the star's cen- 
ter of mass from every grid cell, leaving just the tidal force, 
which ensures the star remains centered in our computational 
domain. Self-gravity is calculated using a multipole expan- 
sion of the star's mass distribution with only the monopole 
and quadrupole moments contributing significantly. 

Different refinement levels have different effective viscosi- 
ties, so it is important that the refinement of cells are chosen 
to match the symmetry of the problem. We use an adaptive 
mesh scheme that refines all zones with p > 10~'^ g cm""* to 
have 8 levels. These zones are refined twice more if they have 
\z\ < 10'" cm, resulting in a total refinement of 10 levels near 
the orbital plane where the shock is expected to form. Each 
block is then divided into 8^ grid cells, making our smallest 
cells 2x10^ cm in width. At pericenter, the star is resolved 
by 10** grid cells. 

Because th e majority of the stars in the loss cone are on 
radial orbits dMagorrian & Tremaind [19991) . the orbit is as- 
sumed to be parabolic. We begin the simulation at f = -10'* 
seconds, the star crosses the tidal radius at f = -10^, and 
reaches pericenter at f = 0. At pericenter, the star is brought 
to within rp = 10'^ cm of the black hole, which corresponds 
to an impact parameter (3 = rt/rp = 7. Orbital energy is lost 
to heat and rotation injected into the star during the passage, 
so we expect some deviation from a parabola. However, for 
large (3 the fractional loss of energy is small. At pericenter, the 
kinetic energy of the star's bulk motion relative to the black 
hole at is PCMqAJ^/ = 2 x lO^-' ergs, while our numerical 
results show only 2 x 10^" ergs are injected as heat and inter- 
nal motions. Although general relativistic effects do affect the 
shape of the orbit somewhat for /3 = 7, we restrict ourselves to 
a purely parabolic encounter in Newtonian gravity. 

3. HYDRODYNAMICS OF STELLAR DISRUPTION 
3.1. Initial Approach 
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Fig. 2. — The hydrodynamics of shock formation in the core of a tidally 
compressed star. The snapshots show the progression of gas temperature dur- 
ing the bounce phase of the passage in a slice that passes through the center 
of mass, is perpendicular to the orbital plane, and parallel to the minor axis of 
the star. The orange coloring in these figures indicates logjQ T where p > 1 g 
cm"' . Regions of subsonic/supersonic collapse/expansion are indicated with 
roman numerals, with solid contours indicating the transition from collapse 
to expansion. The dashed contour is where M- = v-/cs = 1, the mach 1 sur- 
face relative to the z direction. The angle of the minor axis to the x-z plane is 
126.2° in (a), 127.3° in (b), 127.8° in (c), and 141.9° in (d). 

A solar-type star on a loss-cone orbit is initially in hydro- 
static equilibrium as it advances on a parabolic trajectory to- 
wards the black hole, eventually nearing r^, the distance at 
which the black hole's tidal forces are comparable to the star's 
self-gravity. In this region, the velocity and gravitational force 
vectors v and Fg of the star relative to the black hole are al- 
most parallel, resulting in the star being stretched into a pro- 
late spheroid along its direction of travel (Figure[T]). 

As the star continues in its orbit beyond rt, changes angle 
with respect to v, and the star is no longer stretched along its 
major axis. The angle of the star's major axis at is partially 
preserved during the passage. As a result, the star is not per- 
fectly aligned with the trajectory of the center of mass. For 
(3^1, the motion of the fluid after crossing rt (until reaching 
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Tp) is well-described by the trajectories of a system of col- 
lisionless particles. During this phase, the star is vertically 
compressed by the tidal field, leading to a velocity field that 
is directed towards the orbital plane and whose magnitude is 
related to z, the distance above the plane. 

3.2. Rebound 

In the collisionless approximation, every particle would be 
expected to cross through the orbital plane shortly after the 
star's center of mass crossed pericenter However, as the pres- 
sure of the gas in the orbital plane increases, it eventually be- 
comes large enough to overcome the tidal field, which is de- 
creasing as the star moves away from the black hole. When 
the vertical pressure gradient of the gas at the leading edge 
of the star becomes larger than the tidal field, the flow re- 
verses and compression waves propagate outwards in all di- 
rections. The compression waves perpendicular to the orbital 
plane travel through a decreasing density gradient and even- 
tually steepen into shock waves. Shock waves that form on 
the edges of the star are able to reach the surface quickly as 
the star has little vertical extent in these regions. The first of 
these shocks form in the leading edge of the star ^^30 seconds 
after pericenter. Subsequent shock formation sweeps across 
the star (parallel to the star's major axis) at ^ Vp, with the 
edges of the star rebounding before the interior, as illustrated 
in Figure (|2|i. 

The largest densities and pressures are reached when the 
compression wave crosses what was the core of the star prior 
to the encounter. We find that Pmax = 1-3 x lO'^ dyn cm"^ 
at a density p^ax = 310 g cm"^*, occurring 53.1 s after peri- 
center (Figure |2] panel c). These conditions persist for only 
~10 seconds, and substantial nuclear burning can only be trig- 
gered on this timescale for T ^3x 10^ K and p ~ lO'' g cm"-' 
^Champagne & Wiescher 19 921). Compar ison with the results 
presented in IBrassart & Lumineli (l2008h shows that this is 
substantially lower than the predictions of one-dimensional 
simulations. This is expected because the density at the lead- 
ing edge of the star decreases after infall reversal, allowing 
the pressure build-up at the interior to be relieved by leaking 
into the newly rarified post-bounce material. The reduction in 
pressure due to this effect moderates the compression waves 
that form downstream in denser parts of the star by decreas- 
ing the vertical pressure gradient d^P{x,y,t), a feature that is 
absent from one-dimensional simulations, which assume that 
all parts of the star collapse at the same time and ignore the 
effects of neighboring regions. 

Additionally, a consequence of simulating tidal disruptions 
in three-dimensions is degraded linear resolution, which re- 
duces the mid-plane pressure by a factor x ^ where / is 
the width of a grid cell. Applying Brassart and Luminet's 
scaling laws to our progenitor model, we find that they predict 
a compressed core pressure that is x = 100 times larger than 
what is found in our simulation. Because the pressure gradi- 
ent and resolution effects both lead to a decrease in the mid- 
plane pressure, we expect that the actual value of the pressure 
at the mid-plane should fall somewhere between Brassart and 
Luminet's results and our own. 

Shortly after the densest part of the star reaches maximum 
compression, the trailing edge enters its rebound phase. As 
the infalling material is halted by internal pressure at the in- 
terior, the high pressure wave continues to sweep through the 
trailing half of the star at velocity ^ v^. Because the compres- 
sion wave is traveling through a decreasing density gradient, 
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Fig. 3. — Histogram of the energies summed over all fluid elements during 
maximum compression, £, = Kj - Ui where Ki = ^niiCvbulk + vl)^ and Ui = 
GMi,mi/r. Vbulk is the star's bulk velocity, while m; and vi are the rth fluid 
element's mass and velocity. Energies are distributed among 30 bins and the 
data is smoothed using a spline fit. Curves are about 1.25s apart in time, with 
colored arrows indicating effect of the pressure wave propagating through the 
star. 

the lighter material piles up on the denser material towards the 
core, resulting in an increase of velocity of the trailing edge 
relative to the star's center of mass. This has the effect of 
binding more material to the black hole (Figure[3]) because the 
force is opposite to the direction of the star's motion. As a re- 
sult about 5% of the mass that was already bound to the black 
hole becomes more deeply bound after the maximum com- 
pression wave has crossed the entire star One-dimensional 
models do not replicate this effect mainly because they as- 
sume each column of gas collapses independently, neglecting 
both compression and rarefaction waves that propagate from 
other regions of the star. The relative speeds between any two 
parts of the star is proportional to the total spread in veloc- 
ity Av cx Vp cx (lReeslll988b . Because the temperature of 
the st ar is largely uniform and the sound speed is proportional 
to P jLuminet & CarteJl 19861) . a given fluid element will be 
causally connected with a larger fraction of the star for larger 
penetration factors, and thus multidimensional effects become 
increasingly important. 

3.3. Free Expansion 

After the entire star has rebounded, it proceeds to expand 
and adiabatically degrade its internal energy content. Once 
again, the star's pressure becomes unimportant in determin- 
ing its dynamics; as the star leaves the vicinity of the black 
hole, the effects of tidal stretching begin to dominate. As the 
leading edge of the star is slightly closer than it would be on a 
purely parabolic trajectory, the star experiences a torque and 
slowly rotates counter-clockwise when viewed from above 
(Figure IHi. 

Because different parts of the star lie at different dis- 
tances to the black hole, they experience slightly different 
torques and acquire a range of angular velocities. This gen- 
erates a spiral feature similar to that seen in Figure 2 of 
lEvans & Kochanek (1989). This differential rotation leads to 
the star folding onto itself, and as a consequence high pressure 
regions in the interior are lifted towards the star's atmosphere 
and result in the ejection of material (third panel of FigurelU. 
These ejections becomes prominent ^200 seconds after peri- 
center, with the exposed surface area being determined by the 
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exact structure of the differential rotation and the state of the 
remnant's atmosphere. 

Eventually, the star leaves the immediate vicinity of the 
black hole and enters a free expansion phase. As approxi- 
mately half of the star's mass has a negative binding energy 
relative to the black hol e, the bound material w ill eventually 
return to pericenter ( Ramirez-Ruiz & Rosswod 12009), form 
an accretion disk, and feed an AGN phase for a few months 
(lReeslll988h . As we will describe in the next section, the re- 
bound of the remnant also produces a clear observable signa- 
ture, which can be used to predict when a quiescent black hole 
will become active. 



4. SHOCK BREAKOUT 

The immediate signatures of tidal disruption will be observ- 
able in the X-rays, which correspond to a surface temperature 
of ~ 10^ K. Because the ratio between the star's vertical ex- 
tent and its width in the orbital plane is less than 1:10, the 
light curve produced by shocks propagating perpendicular to 
the orbital plane should have roughly a cos(V-') dependence, 
where -if) is the angle between the orientation of the orbital 
plane and the line of sight. Due to the complex nature of the 
shocks generated by this event, the escaping radiation does 
not necessary lead to a single burst. In fact, our simulation 
shows two distinct peaks (Figure|5]i. 

Because our simulation only resolves the star out to den- 
sities of ~ lO"-' g cm""*, our calculated surface lies beneath 
the true t = 1 surface for Thomson scattering in a stellar at- 
mosphere. This surface corresponds to densities of nip/aoil = 
10"^ g cm""* in our simulation, where / is the length of the 
smallest grid cells. Additionally, the simulation treats the en- 
tire star as a = 5/3 ideal gas, which is clearly not true in the 
tenuous outer layers of the star's atmosphere that are radiation 
pressure dominated for temperatures of 10^ K and densities 
of P ^ 1 g cm"-'. Because of these limitations, we must ex- 
trapolate from the conditions at the surface of the simulation 
to determine the properties of the post-shocked atmosphere. 

We define the surface of our simulation to be where Cgas > 
5ei-ad such that the gas temperature is still representative of 
the energy content of a given fluid element. As the pho- 
tons below the r = 1 surface are trapped on the shock cross- 
ing timescale, we can assume the shock has no radiative 
losses and the properties of the post-shock material are ac- 
curately de scribed by the Rankine-Hugoniot jump conditions 
( lShulll992h . For a perfect gas of ionized Hydrogen the en- 
ergy density is Cgas = 3nki,T, where n is the number density 
of atoms and T is the pre-shock gas temperature. Thus, the 



post-shocked material will have e : 



■ky,T /nipj, where 



Vsh is the shock velocity. This expression is controlled by the 
v^j, term in the star's outer layers where T is small, and thus 
the energy density of a post-shocked region is oc pv^j^, the ram 
pressure of the post-shock material. 

If we take the energy density ej = 3>pJibTsm~^ (where the 
subscript s refers to the surface conditions) to be equal to the 
post-shock energy density at the base of the atmosphere, we 
can then calculate the energy density e as a function of height 
z above this surface by considering the self-similar solutions 
of a shock propagating through a decreasing density gradi- 
ent. In such solutions, the evolution of the shock velocity as it 
moves towards the photosphere is determined by 7ad and the 
dependence of p on z, where z is the distance from the base of 
the atmosphere. 



4.1. Self-Similar Solutions 

Solutions for the post- shock conditions in the stron g shock 
limit of both expo nential (^ Grover & H ardy 1966: Haves 19681: 
IZel'Dovich & Raizeri. 19671) and power-law (Sakurai I960) at- 
mospheres are readily available. The outer layers of an F = 3 
polytrope are well- described by a power-law distribution with 
index 1.5 < 5 < 3 (iMatzner & McKedlI999h : we assume that 
(5=1.5 for the remainder of this work. For a power-law atmo- 
sphere where z is the distance from the base of the atmosphere 
and liQ is the pressure scale height, the density p{z) is 



p(z) = pA^+zK^y 



(2) 



If I'sh.s = Ph is the shock velocity at the base of the atmo- 
sphere, then v(z) is 



v(z)=^[l+z/j-'](i-)/° 

7ad+l 



(3) 



Because the post-shock energy density is just equal to the ram 
pressure in the strong shock limit, the total energy density e(z) 
is simply 

e(7) = es[l+z/jo']2"-"'/"-^ (4) 

At the base of the atmosphere thermal pressure dominates, 
but as the density of the gas decreases with z more rapidly than 
the temperature, the atmosphere becomes radiation pressure 
dominated within a couple scale heights and thus 7ad = 4/3. 
For our assumed values of 8 and 7ad, Sakurai's method leads 
to a similarity exponent a = 0.7774. 

As the atmosphere is strongly compressed by the time- 
integrated tidal gravitational field, its vertical scale at shock 
breakout is significantly reduced. To estimate the atmo- 
sphere's size at breakout, we assume that collapse is self- 
similar. Because most of the stretching in a close encounter 
takes place before the star reaches rp, the cross-sectional area 
A* should not be a strong function of /?, provided that P ^ I. 
Thus, the star's vertical size is directly related to the change 
in the star's volume 

^ = ^ (5) 

where H.; (= and Hf are the initial and final heights of the 
star above the orbital plane and Vi and Vf are the initial and 
final stellar volumes. The vertical velocity vj^ of the stellar 
material just prior to rebound is approximately equal to /3cs 
((Carter & Luminet 1983), and because nearly all of the in- 
ternal motions during free-fall are vertical, the total kinetic 
energy of the star can be estimated as 
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GMl 



(6) 



At rebound, the kinetic energy of the infalling material is 
converted into internal energy. Because the star's original in- 
ternal energy is U, ~ GMI/R^,, the initial and final internal 
energies are simply related, Uf = 0^Ui. Assuming the com- 
pression is adiabatic. 



V o^U '-"'ad . 



and thus 



Hf = Hi(| 



= Hi/3' 



(7) 



(8) 
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Fig. 4. — Three snapshots showing a sHce of the star in the z = plane, where t is time after pericenter. The white arrow in each panel indicates the direction to 
the black hole's center The orange coloring of the panels indicate logigP where p > 10^^, where P,ni„ and Pmax are set to the minimum and maximum pressure 
in each snapshot. In panels (b) and (c), contours of constant angular velocity lo = v„-/\r\ are drawn in cyan showing deviations from ligid-body rotation. 



For 7a(j = 5/3, equation dHJ gives Hf = TiiP ^. In a self- 
similar collapse, the scale height /ib will be reduced by a factor 
Hf/Hi, and thus /ib = IiqP'^ . 

Because the velocity of a given layer is position dependent, 
Av is non-zero and the post-shock energy density e is a de- 
creasing function of t. The energy output will therefore be 
controlled by a balance between photon diffusion and adia- 
batic expansion, with most of the photons being released from 
a layer where these two timescales are comparable. Assuming 
plane-parallel geometry, these timescales can be estimated for 
a layer at a given z as 



Tad ■ 



Tdiff : 



erfl/V2 -1 



dv(z) 



t(z) p(z)(Tlh 



(9) 
(10) 



where CTth is the Thomson cross-section and nip is the proton 
mass. Equating these two expressions and solving numeri- 
cally for z determines Zmin, the deepest layer in the atmosphere 
that will contribute significantly to the emission. 

With the compressed density profile, shock jump condi- 
tions, and depth of the emitting region determined, we can 
calculate the luminosity and spectrum of shock breakout. As 
the energy density of the post-shocked material is primarily 
given by the energy density of radiation, we can estimate the 
photon temperature in the atmosphere as T^^ = e(z) / a^, where 
Ub is the Stefan-Boltzmann constant. While the temperature 
does decrease as z approaches H (the distance from the base 
of the atmosphere to the photosphere), this decrease is small 
because T cx e'/"*. 

The spectrum of photons near the t = 1 surface may be 
significantly non-thermal because the photon spectrum takes 
a non-ne glible amount of time to approach a Planckian dis- 
tribution ^KatFetaTIIJOOll). Additionally, the shock velocity 
near the surface is ^ 0.03c, and special relativistic effects may 
also contribute to a non-thermal spectrum. While these effects 
would le ad to the producti on of harder photons that are easier 
to detect ( iBand et al.l2008h . we conservatively assume that the 
observed photon spectrum is a combination of blackbodies. 

The total energy released by the event can be estimated by 
integrating the energy density over z and multiplying by twice 



the cross-sectional of the star at pericenter 2A, = lirab: 



Etot = lirab 



e(z)dz, 



(11) 



where a and b are the semi-major and semi-minor axes of the 
star We can define an effective height parameter Z 



Z = 



2(1-q)/q-i5 



- h Q; / _z_ 



dz 



(12) 



which is independent of position on the surface of the star if 
we assume that the atmosphere is always perpendicular to the 
orbital plane and that h\, and S are constants. The average 
luminosity of the event is then given by substituting Z into 
equation JTTJ and dividing by the length of time for the star's 
center of mass to cross pericenter 



■ 2a/vp 



6R 



3/2 



1/3 ,,1/6' 



which yields 



L = TTVp^JEsZ, 



(13) 



(14) 



where we used equation (|4|l to substitute for e(z) and our defi- 
nition for Z. Evaluating equation (fT4l) for the average found 
at the surface of our simulation, assuming the initial scale 
height has the solar value of ho ~ 1.4 x 10^ cm, and using 
the measured crossing time ?„- = 30 s, we find Z = 3 x lO'*'* 
ergs s"' . 

While our estimation determines the average luminosity Z, 
calculating the luminosity as a function of time is more dif- 
ficult because is a function of position on the surface. Be- 
cause a detailed record of the surface conditions are produced 
by our simulation, we can calculate what the luminosity and 
spectrum will look like as a function of time. Most of the lu- 
minosity is released from the deepest contributing layer of the 
atmosphere e(Zmin), where the diffusion time to the surface is 
~ 10"^ s. This timescale is significantly shorter than so we 
approximate each grid cell's contribution to the luminosity as 
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Fig . 5 . — Light curve produced by the breakout of shocks across the surface 
of the disrupted star, where s is the flux in photons per second and t is time 
since the star's center of mass crossed pericenter Photons are binned into 
categories based on energy. The dotted hues show the raw output of the 
luminosity calculation for a given bin, while the solid Unes show a 1 s moving 
average of the same data. 



a delta function 



L{x,y,t) = 



At 



(15) 



where / is the grid cell size in cm. At is the time increment 
between two data dumps, and t\,{x,y) is the breakout time for 
a given location on the surface. We use At = 0.1s, but the 
value chosen only affects the smoothness of the resultant light 
curve so long as Af ^ ta- Because the shock front is not al- 
ways easy to detect in our simulation due to limited resolution 
in the perpendicular direction, t\, is set by when a given area 
element's energy density peaks. The luminosity as a function 
of time is then simply a sum of the individual contributions of 
the surface grid cells 



(16) 



where e^.i and fb,, are the energy density and breakout time of 
a given cell, respectively. 

Equation ( fT6l ) allows us to generate a light-curve for our 
disruption event (Figure|5]). We assume that each grid cell pro- 
duces a blackbody spectrum of photons defined by the photon 

1 /4 _i 

temperature of that cell, Tph,,- = ejj . By summing these 
spectra, we can obtain a photon distribution at time t. Note 
that there are two peaks at ^^30 and ^^50 seconds correspond- 
ing to shocks that originate from the limbs of the star. As 
discussed in |3] this double peak feature arises because com- 
pression waves that originate near the star's center of mass 
tend to travel parallel to the orbital plane in the direction of 
the star's orbital motion into regions that have akeady begun 
expanding. This compression moderation results in the non- 
production of shocks in this region, and thus there is a period 
of relatively low luminosity as the shocks cross the star's cen- 
ter of mass. The time difference between the two maxima in 
the light curve is approximately equal to fcr (eq- Ell. 

Our post-shock surface conditions show peak breakout tem- 
peratures in excess of 8 x 10^ K at densities of 30 g cm"-'. 
iKobavashi et al.. (.2004) estimated that the /3 = 5 event will 



yield photons of energy of 2.2 keV on average, which is in 
agreement with our mean surface temperature. 

5. OBSERVABILITY 

5 . L X-ray Transient 

The duration, color, and luminosity of a disruption breakout 
event depends on /3, Mh, and the parameters of the star being 
disrupted. Ideally, one would want to construct a predictive 
model for TSBs by performing a detailed simulation for each 
possible combination of parameters. However, a thorough ex- 
ploration of the full parameter space would require many sim- 
ulations similar in scope to the simulation presented in this 
work. Fortunately, progress can still be made as our simula- 
tion provides an accurate benchmark for the peak luminosity 
Lpeak of a deeply-penetrating breakout event. The results of 
the simulation can then be used together with some simple 
scaling arguments to construct a function that can describe the 
features of disruption breakouts for a variety of encounters. 

We can gain some insight by considering the scaling of the 
various parameters that the luminosity (eq. fl4l i depends on. 
The pericenter velocity is dependent on the black hole mass 
and the closest-approach distance, which gives a scaling of 

Vp (X (3^I^m\^^ . For = 5/3 (appropriate for the mid-plane), 
pec 0^ and T cx (iLuminet & CarteJll986h , which implies 
that 65 cx pv^ oc (3^ . Because the total height of the atmosphere 
H scales as /?"■', large /3 events are reduced in output by the 
decrease in volume of the emitting region. For ultra-close 
encounters the escape velocity can be comparable to c, and 
thus the emitted photons are gravitationally redshifted, with 
the energy of each being divided by a factor of (l + ^r^/r^Y 
Additionally, the passage of time in the star's rest frame is 
slower than that of an observer for which r ^ rj by a fac- 
tor of ^1 -Ts/rp, leading to a further decrease in luminosity. 
Accounting for these effects, the peak luminosity of a TSB 
roughly scales as 



-^peak * 



1/3 y/l-^Ap 
2 + rs/rp 



(17) 



A more careful calculation considering the detailed proper- 
ties of the shocked, self-similar atmosphere reveals that L also 
depends on the choice of 5. The dependence on (3 and Mh still 
roughly follows the scaling of equation ( fTTI i. but because the 
volume of the emitting region depends on where r^d = Tjad 
(eqs.|9]and[T0li, the luminosity of a given event must be com- 
puted numerically. The full numerical solution for Lpeak of a 
(5 = 1.5 power-law atmosphere over a range of j3 and Mh is 
shown in the left panel of Figure 

We now follow the procedure of TWang & MerrittI (120041) to 
estimate the rate of detectable events. Because the tidal ra- 

1 /3 

dius grows as and the Schwarzschild radius grows as Mj,, 
disruptions with /3 > 3 are only possible for black holes with 
^ IO'^Mq. This means a significant detection rate is only 
obtained when considering the l ow end of the galaxy lumin os 
ity function, as characterized bv lTrentham & Tullvl ( 120021) 



_jq-0.4(A1-A1je) 



, CtdE+l 



(18) 



X e 



dM. 



In this expression, A'dE.o is a normalized dwarf elUptical 
galaxy number density, A^dE is the cutoff magnitude for the 
Schechter function, and OfdE describes the faint-end slope. As 
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Fig. 6. — Peak luminosity and detectability of events as a function of Mj, 
and (3, assuming that the outer atmosphere is described by a 5 = 1.5 power 
law. The left panel shows logigLpe^k for a given event, while the right panel 
shows the number of expected detections by EXIST fGrindlay 2004) , nor- 
malized to the maximum rate. Our simulation parameters are indicated by 
the white circles, while the black crosses indicate the location of the events 
with the largest predicted detection rates. Disruptions are constrained by the 
/3 = /3,Tiax curve (shown as a dashed line), which indicates the largest /3 for a 
given black hole mass where rp > rs . 



in I Wang & Merritj ( 120041) . we want to write equation ([TS]! as a 
functio n of Mh. We use the scaHng relation of Magorri an et alj 
(119981) to write the absolute magnitude A4 in terms of bulge 
mass Mhiiire, which is sim ply related to Mh = 10~^ '"Mbuige 
dMerritt & Ferrares^l2001h . We assume that only nucleated 
dwarf elliptical (dEn) galaxies contain black holes, with the 
nucleated fraction F„ scaling linearly with Ai. To properly 
scale A^dE.o, we assume that the space-averaged dE density is 
equal to the number of dEs in the Virgo cluster spread into a 
sphere with radius equal to the distance of the Virgo cluster 
(approx. 0.2 Mpc"-'). 

Because we are in the regime where rp < C r,, we know that 
the disruption rate scales linearly with rp (lReesl[T988h . and 
Wang's expression for the disruption rate becomes 



N = 6.5x IQ-'* 



yr 



lO^Mm 



-0.25 



/3- 



(19) 



The maximum distance Rd to a detectable event can be ex- 
pressed as a function of instrumental sensitivity as 



ttL, 



peak 



hv 



(20) 



in which Fj{v) is the burst sensitivity as a function of fre- 
quency, the TSB spectrum is given by the Planck function 
at a temperature Tph oc fi^^'^, h is Planck's constant, <Tb is the 
Stefan-Boltzmann constant, and v^i^ and z^^ax are the mini- 
mum and maximum frequencies accessible to the instrument. 
An instrument is sensitive to all events contained within the 
volume 

^ (21) 



We now have all the pieces needed to determine a detection 
rate in terms of /3 and Mh. By multiplying the number density 
of black holes NaeFu (eq- [TST l by the disruption rate .^V for a 
black hole of a given mass and penetration factor (eq. [T9l l, 
and by the volume Vd in which those events are detectable 
(eq. l2lT i. we obtain the following expression for the total rate 



• E + dE + Bulges 
i AGN 

• Globular Clusters 
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Velocity Dispersion a (km s ') 

Fig. 7. — Bulge velocity dispersion (Tc vs. BH masses in E, dE, and bulge- 
dominated spiral galaxies (circles), AGN (triangles), and globular clusters 
(pentagons). The colored contours represent the detection rate relative to the 
maximum detection rate as in Fig. j6), averaged over all (3, with the max- 
imum indicated by the dashed line. Plot mo dified from Geha et al. 1 2()03) , 
including the BH masses compiled bv lNovola et al.i 120081) . 



of detection 



47r 



NNdEFnVidMudP, 



(22) 



Mm 



in which f2 is the solid angle covered by the instrument's 
field of view and M,„ax = {c-RQ/AGMQ/Sf'^ 



at which rp < r^, the Schwarzschild radius. 



is the mass 
The integral 



over P terminates at /3 = 3 because shocks have not been 
seen in one-dimensional ca lculations for passages with /3 < 3 
(iBrassart & Lumine3l2008h . However, our simulation shows 
that shocks form in the limbs of the star where conditions are 
substantially different than the core, and thus shocks may still 
be produced for < 3 in some cases. 

The detection rate is evaluated using the parameters of the 
proposed EXIST telescope, which has spectral coverage from 
3 to 1000 keV, an average burst sensitivity of Fj = 0.2 cm~^ 
s"', a nd - 20% sky coverage dGrindlavl l2004t iBand et alj 
120081) . Using these values, equation ( |22] ) predicts that EX- 
IST should detect approximately 1 TSB per year We can also 
use the argument of integrals of this equation to generate a 
map in (/3,Mh) space (Figure|6] right panel) to predict which 
events will generate the most detections. Comparison with 
the left panel of the same figure shows that these events have 
a characteristic luminosity of ~ 5 x 10'*'* ergs. Despite being 
similar to the Eddington value of lO'^'^Mf, for a SMBH, detec- 
tion is infrequent because of the short duration of the shock 
breakout. By integrating equation (l22l i over all (3, we can also 
determine how black holes of different masses contribute to 
the detection rate (Figurel?), with the peak rate corresponding 
to a black hole mass of 2 x IO^Mq. The predicted peak in the 
detection rate is predicated on the assumption that SMBHs 
obey the black hole to bulge mass relation for Mh < IO^Mq. 
Therefore, the detection of TSBs would test the validity of 
this assumption and potentially provide compelling evidence 
for SMBHs in low mass galaxies. 

5.2. Gravitational Wave Signal 
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In addition to the X-ray breakout signature, gravitational 
waves are also expected to radiate from the encounter, mostly 
originating from the changing location of the star relative to 
the black hole. We can estimate the strength of these waves 
by approximating the second time derivative of the moment 
of inertia tensor as ~ GMR^ / c^P^ . If we assume that the 
period P ^ r^/vj, = r^/ GM^^, the length scale of variation R ~ 
Tp, and the mass M = M*, we find that the gravitational wave 
ampUtude h is (iKobavashi et al.ll2004 



/3 



(23) 



in which d is the distance to the event. For our (3 = 7 simula- 
tion, we calculate that li ^ 10"^' for li = 10 Mpc. Gravitational 
waves will also be radiated as the star itself changes shape, but 
these distortions are far smaller than those generating by the 
changing quadrupolar moment of the star-SMBH system. As- 
suming that R ^ the gravitational wave amplitude is 



P'G^Mj 
dc*R^ 



(24) 



Note that this expression is independent of Mj,. For d = 10 
Mpc, Ji is 10"^^ for our encounter, beyond LISA'S sensitivity 
range. Because this expression is proportional to P^, the com- 
pression of the star can be a substantial contribution to the 
gravitational wave signature for very deep passages (/3 > 25). 

If both a gravitational wave signal and a TSB signal are 
available for the same event, additional information about the 
encounter can be obtained. While both observational signa- 
tures each place upper and lower limits on the properties of 
the disrupted star, either signature on its own cannot uniquely 
constrain the star's characteristics. By assuming a stellar 
mass-radius relationship M^, = MqR^, equation i23[ can be 
combined with equation (fT3T l to calculate the mass of the dis- 
rupted star 



(25) 



in which can be estimated by measuring the distance be- 
tween the two peaks in the light curve (e.g. Figure |5]l. By 
determining the masses of disrupted stars, a distribution of 
stars that occupy the loss cone can be derived. As the IMF 
in the vicinity of SMBHs is poorly characterized even in our 




own galaxy (lAlexandeiil2005[) . a coincident detection of both 
the gravitational and TSB signatures would, for the first time, 
allow us to investigate the IMF in close proximity to extra- 
galactic SMBHs. 

6. CONCLUSION 

While the simulation presented in this work features the 
highest resolution of the tidal disruption of star in 3D to 
date, it still has a few shortcomings. Our model does not 
have enough linear resolution to resolve the sharp pressure 
gradients that develop in deeply penetrating encounters, and 
thus the mid-plane pressure in the rebound phase is certainly 
underestimated by some factor In addition, for close pas- 
sages such as ours, GR effects start to become important. For 
5 < /3 < 10, the orbit is better characterized by a Paczynski- 
Wiita potential, but for /? > 10, a ful ly general relativis tic 
treatment of the Schwarzschild metric ( iFrolov et alJI 19941) is 
required. The orbits are then not ellipses, but may have two 
or more pericenter transversals and, as a result, could lead to 
the formation of multiple shocks. A full understanding of the 
compression process in such cases will require detailed GR 
hydrodynamical simulations. 

In this work, we have provided numerical details of how 
a solar-type star is stretched, squeezed, and strongly shocked 
during an encounter with a massive black hole. We then cal- 
culated the radiation a distant observer might detect as the 
observational signature of the accompanying shock breakout. 
If detected, an L ~ 10"*"* erg/s burst that fades within a few 
minutes and exhibits the predicted double-peaked signature 
in the soft X-rays would be compelling testimony that a star 
experiencing an ultra-close encounter with a black hole can be 
disrupted and compressed to such an extent that shock waves 
can be triggered. 
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APPENDIX 



TABLE 1 
Glossary of symbols. 



Term 


Meaning 


Term 


Meaning 


Fundamental constants 




Shock breakout 






Radiation constant 


ho 


Initial atmospheric scale height 


c 


Speed of light 


K 


Atmospheric scale height at shock-breakout 


G 


Newton's constant 


H 


Depth of compressed atmosphere 


h 


Planck's constant 


n 


Atomic density 


kb 


Boltzmann's constant 


% 


Shock-breakout time 


nip 


Proton mass 


T 


Gas temperature 


Mq 


Solar mass 




Photon temperature 


Rq 


Solar radius 


V 


Post-shock velocity of atmosphere 




Stefan-Boltzmann constant 




Shock velocity 


^th 


Thomson cross-section 




Depth of deepest emitting layer 






Z 


Effective height 


Encounter parameters 




a 


Similarity exponent 


a 


Semi-major axis of elongated star 


7!id 


Adiabatic index 


A, 


Cross-sectional area of elongated star 


5 


Power-law index of stellar atmosphere 


b 


Semi-minor axis of elongated star 


€ 


Total energy density 


H 


Star's extent above orbital plane 




Ideal gas energy density 


Mh 


Black hole mass 


^rad 


Radiation energy density 


Me 


Black hole mass in 10^ Mq 


p 


Mass density 


M, 


Mass of star 


T 


Optical depth 


Mo 


Stellar mass-radius normalization 


fad 


Atmosphere expansion timescale 


''s 


Schwarzschild radius 




Radiation diffusion timescale 


n 


Tidal radius 






fp 


Distance at pericenter 


SMBH distribution 




R, 


Initial stellar radius 




Nucleated fraction of dE galaxies 


/ 


Time relative to time of pericenter 


^bulge 


dE bulge mass 


fcr 


Time for star to cross pericenter 


M 


Galaxy magnitude 


u 


Total internal energy of star 


■MdE 


dE Schechter function cutoff magnitude 


H 


Velocity of star at pericenter 


NdE 


dE number density 




Velocity of star's vertical collapse 




dE Schechter function faint-end slope 


V 


Volume of star 






/3 


Penetration factor 


Observability 




r 


Polytropic index of star 




Planck distribution 




Stellar mass-radius relation exponent 


d 


Distance of disruption from observer 






h 


Gravitational wave ampHtude 


Simulation parameters 




Ft 


Burst sensitivity of instrument (ph s"' cm"^) 


I 


Width of smallest grid cells 


L 


Bolometric luminosity 




Pressure maximum in simulation volume 


L 


Average bolometric luminosity 


T, 


Gas temperature at simulation surface 




Peak bolometric luminosity 


At 


Time between simulation plot files 


N 


Tidal disruption rate of a given SMBH 


Vs 


Post-shock velocity at simulation surface 


Rd 


Maximum detectable distance 


Vsh,s 


Shock velocity at simulation surface 


Vd 


Volume of sphere with radius i?j 




Distance from orbital plane 


n 


Solid angle coverage of instrument 




Total energy density at simulation surface 






Pmax 


Mass density maximum in simulation volume 






Ps 


Mass density at simulation surface 






X 


Pressure reduction due to hmited resolution 







